function plot_FOX_eta_m_bs_vMixed(Fox_model,subsample,inc,set_seed,matlab_seed,get_nb_draws)

% Collect_Fox_discrete_FEdemoshrt_grID_wtp_tomlab_varagin_j_bs(subsample,inc,set_seed,matlab_seed,get_nb_draws)
% do \\c3\rdat\SHoude\Research\EEgap\EEgap_scripts\plot_Fox_contour_3dim_bs.do 

path_to_results_folder='/Users/shoude/Dropbox/eegap/EEgap_data_code_heter_SM/Results_Estimation/';
 
param_results_fox = strcat('plot_jointpdf_Fox_',Fox_model,'_mixedFEdemoshrt_WTP_v',num2str(subsample),'_inc',num2str(inc),'_sd',num2str(set_seed),'_cty_tau1_grid_v3_all_bs.csv');
param_est     	  = csvread([ path_to_results_folder param_results_fox]);
	
p10 = [1:10]/10;
pdf10_edge = quantile(param_est(:,3),p10);
pdf3_edge =[0 5 15 max(param_est(:,3))];
[pdf3_color_id] = discretize(param_est(:,3),pdf3_edge);

tstats2_edge = [0 1.1 max(param_est(:,7))];
tstats2_color_id = discretize(param_est(:,7),tstats2_edge);
param_est_g=[param_est, pdf3_color_id, tstats2_color_id];

group_pdf=findgroups(pdf3_color_id,tstats2_color_id);
nb_group=max(group_pdf);
param_est_g=[param_est_g, group_pdf];
for j=1:nb_group
	group_id(:,j)= group_pdf==j;
end
%pdf/tstats/group
% 1 1 1
% 1 2 2
% 2 1 3
% 2 2 4
% 3 1 5
% 3 2 6
 fig1=figure(1)
hold on
% pdf < 5, tstats < 1.1
h1=scatter3(param_est(group_id(:,1),1),param_est(group_id(:,1),2),param_est(group_id(:,1),3),50,'o');
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0.5];

% pdf < 5, tstats >= 1.1
h2=scatter3(param_est(group_id(:,2),1),param_est(group_id(:,2),2),param_est(group_id(:,2),3),50,'*');
h2.MarkerFaceColor = [0 0 0];
h2.MarkerEdgeColor = [0 0 0.5];

% pdf [5,15], tstats < 1.1
h3=scatter3(param_est(group_id(:,3),1),param_est(group_id(:,3),2),param_est(group_id(:,3),3),100,'o');
h3.MarkerFaceColor = [0 0 0];
h3.MarkerEdgeColor = [0 0.75 0];

% pdf [5,15], tstats > 1.1
h4=scatter3(param_est(group_id(:,4),1),param_est(group_id(:,4),2),param_est(group_id(:,4),3),100,'*');
h4.MarkerFaceColor = [0 0 0];
h4.MarkerEdgeColor = [0 0.75 0];

% pdf > 15, tstats < 1.1
h5=scatter3(param_est(group_id(:,5),1),param_est(group_id(:,5),2),param_est(group_id(:,5),3),150,'o');
h5.MarkerFaceColor = [0 0 0];
h5.MarkerEdgeColor = [1 0 0];

if inc ~= 1
% pdf > 15, tstats >= 1.1
h6=scatter3(param_est(group_id(:,6),1),param_est(group_id(:,6),2),param_est(group_id(:,6),3),150,'*');
h6.MarkerFaceColor = [0 0 0];
h6.MarkerEdgeColor = [1 0 0];
end

xlabel('Coefficient on Price: \eta','Fontsize',16)
 ylabel('Misperception: m','Fontsize',16)
 zlabel('Estimated PDF Weight','Fontsize',16)
ylim( [-3,6] );
xlim( [-35,5] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12) 
%a = get(gca,'XTickLabel');
%set(gca,'XTickLabel',a,'FontName','Times','fontsize',18)
 
%Extract mean from FKRB estimator

if subsample==11000
	if inc==16
	 line([-5.431913,-5.431913],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.7650573,0.7650573],[0,0] ,'Color','k' );
	elseif inc==1
	 line([-6.421407,-6.421407],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.3539926,0.3539926],[0,0] ,'Color','k' );
	elseif inc==2
	 line([-6.261386,-6.261386],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.4845289,0.4845289],[0,0] ,'Color','k' );
	elseif inc==3
	 line([-5.604053,-5.604053],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.7280716,0.7280716],[0,0] ,'Color','k' );
	elseif inc==4
	 line([-5.082855,-5.082855],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.8185474,0.8185474],[0,0] ,'Color','k' );
	elseif inc==5
	 line([-5.111113,-5.111113],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.9405871,0.9405871],[0,0] ,'Color','k' );
	elseif inc==6
	 line([-4.592757,-4.592757],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[1.000794,1.000794],[0,0] ,'Color','k' );
	end
elseif subsample==44000
	if inc==16
     if strcmp(Fox_model,'discrete')    
	  line([-5.431913,-5.431913],[-3,6],[0,0],'Color','k' );
	  line([-35,5],[0.7650573,0.7650573],[0,0] ,'Color','k' );
     elseif strcmp(Fox_model,'discretev5') 
      line([-5.357534,-5.357534],[-3,6],[0,0],'Color','k' );
	  line([-35,5],[0.7911761,0.7911761],[0,0] ,'Color','k' );
	 elseif strcmp(Fox_model,'discretev6') 
      line([-5.435024,-5.435024],[-3,6],[0,0],'Color','k' );
	  line([-35,5],[0.7624894,0.7624894],[0,0] ,'Color','k' ); 
     elseif strcmp(Fox_model,'discretev7')    
      line([-5.237717,-5.237717],[-3,6],[0,0],'Color','k' );
	  line([-35,5],[0.7336762,0.7336762],[0,0] ,'Color','k' );
     elseif strcmp(Fox_model,'discrete_3tau')
      line([-5.58872,-5.58872],[-3,6],[0,0],'Color','k' );
	  line([-35,5],[0.7778516,0.7778516],[0,0] ,'Color','k' ); 
     end
    elseif inc==1
	 line([-6.421407,-6.421407],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.3539926,0.3539926],[0,0] ,'Color','k' );
	elseif inc==2
	 line([-6.261386,-6.261386],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.4845289,0.4845289],[0,0] ,'Color','k' );
	elseif inc==3
	 line([-5.604053,-5.604053],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.7280716,0.7280716],[0,0] ,'Color','k' );
	elseif inc==4
	 line([-5.082855,-5.082855],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.8185474,0.8185474],[0,0] ,'Color','k' );
	elseif inc==5
	 line([-5.111113,-5.111113],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[0.9405871,0.9405871],[0,0] ,'Color','k' );
	elseif inc==6
	 line([-4.592757,-4.592757],[-3,6],[0,0],'Color','k' );
	 line([-35,5],[1.000794,1.000794],[0,0] ,'Color','k' );
	end
end


if inc==16
  if strcmp(Fox_model,'discrete')   
         Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
                'pdf: [0,5), t-stat: [1.1,\infty)',...
                 'pdf: [5,15), t-stat: [0,1.1)',...
                  'pdf: [5,15), t-stat: [1.1,\infty)',...
                   'pdf: [15,100], t-stat: [0,1.1)',...
                    'pdf: [15,100], t-stat: [1.1,\infty)',...
                    'mean \eta: -5.43',...
                    'mean m: 0.77',...
                    'Location','NorthWest');
        EG = findobj(Lobj,'type','text');
        LEG.FontSize = 10;
  elseif strcmp(Fox_model,'discretev5')
        Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
                'pdf: [0,5), t-stat: [1.1,\infty)',...
                 'pdf: [5,15), t-stat: [0,1.1)',...
                  'pdf: [5,15), t-stat: [1.1,\infty)',...
                   'pdf: [15,100], t-stat: [0,1.1)',...
                    'pdf: [15,100], t-stat: [1.1,\infty)',...
                    'mean \eta: -5.34',...
                    'mean m: 0.79',...
                    'Location','NorthWest');
        EG = findobj(Lobj,'type','text');
        LEG.FontSize = 10;
        
        %-5.357534	0.7911761	0.153497	0.0747381

    elseif strcmp(Fox_model,'discretev6')
        Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
                'pdf: [0,5), t-stat: [1.1,\infty)',...
                 'pdf: [5,15), t-stat: [0,1.1)',...
                  'pdf: [5,15), t-stat: [1.1,\infty)',...
                   'pdf: [15,100], t-stat: [0,1.1)',...
                    'pdf: [15,100], t-stat: [1.1,\infty)',...
                    'mean \eta: -5.44',...
                    'mean m: 0.76',...
                    'Location','NorthWest');
        EG = findobj(Lobj,'type','text');
        LEG.FontSize = 10;      
        
  elseif strcmp(Fox_model,'discretev7')
        Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
                'pdf: [0,5), t-stat: [1.1,\infty)',...
                 'pdf: [5,15), t-stat: [0,1.1)',...
                  'pdf: [5,15), t-stat: [1.1,\infty)',...
                   'pdf: [15,100], t-stat: [0,1.1)',...
                    'pdf: [15,100], t-stat: [1.1,\infty)',...
                    'mean \eta: -5.24',...
                    'mean m: 0.73',...
                    'Location','NorthWest');
        EG = findobj(Lobj,'type','text');
        LEG.FontSize = 10;
      
     %-5.237717	0.7336762	0.1490462	0.0670671
  
  elseif strcmp(Fox_model,'discrete_3tau')
  Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
                'pdf: [0,5), t-stat: [1.1,\infty)',...
                 'pdf: [5,15), t-stat: [0,1.1)',...
                  'pdf: [5,15), t-stat: [1.1,\infty)',...
                   'pdf: [15,100], t-stat: [0,1.1)',...
                    'pdf: [15,100], t-stat: [1.1,\infty)',...
                    'mean \eta: -5.59',...
                    'mean m: 0.78',...
                    'Location','NorthWest');
        EG = findobj(Lobj,'type','text');
        LEG.FontSize = 10;
        
    %-5.58872	0.7778516	0.1796701	0.0983592
      
  end

elseif inc==1
    Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
        'pdf: [0,5), t-stat: [1.1,\infty)',...
         'pdf: [5,15), t-stat: [0,1.1)',...
          'pdf: [5,15), t-stat: [1.1,\infty)',...
           'pdf: [15,100], t-stat: [0,1.1)',...
            'pdf: [15,100], t-stat: [1.1,\infty)',...
            'mean \eta: -6.42',...
            'mean m: 0.35',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==2
 
 Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
        'pdf: [0,5), t-stat: [1.1,\infty)',...
         'pdf: [5,15), t-stat: [0,1.1)',...
          'pdf: [5,15), t-stat: [1.1,\infty)',...
           'pdf: [15,100], t-stat: [0,1.1)',...
            'pdf: [15,100], t-stat: [1.1,\infty)',...
            'mean \eta: -6.26',...
            'mean m: 0.48',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
 
elseif inc==3
 
 Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
        'pdf: [0,5), t-stat: [1.1,\infty)',...
         'pdf: [5,15), t-stat: [0,1.1)',...
          'pdf: [5,15), t-stat: [1.1,\infty)',...
           'pdf: [15,100], t-stat: [0,1.1)',...
            'pdf: [15,100], t-stat: [1.1,\infty)',...
            'mean \eta: -5.60',...
            'mean m: 0.73',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
 
elseif inc==4
 
 Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
        'pdf: [0,5), t-stat: [1.1,\infty)',...
         'pdf: [5,15), t-stat: [0,1.1)',...
          'pdf: [5,15), t-stat: [1.1,\infty)',...
           'pdf: [15,100], t-stat: [0,1.1)',...
            'pdf: [15,100], t-stat: [1.1,\infty)',...
            'mean \eta: -5.08',...
            'mean m: 0.82',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

elseif inc==5
 
 Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
        'pdf: [0,5), t-stat: [1.1,\infty)',...
         'pdf: [5,15), t-stat: [0,1.1)',...
          'pdf: [5,15), t-stat: [1.1,\infty)',...
           'pdf: [15,100], t-stat: [0,1.1)',...
            'pdf: [15,100], t-stat: [1.1,\infty)',...
            'mean \eta: -5.11',...
            'mean m: 0.94',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

elseif inc==6
 
 Lobj=legend('pdf: [0,5), t-stat: [0,1.1)',...
        'pdf: [0,5), t-stat: [1.1,\infty)',...
         'pdf: [5,15), t-stat: [0,1.1)',...
          'pdf: [5,15), t-stat: [1.1,\infty)',...
           'pdf: [15,100], t-stat: [0,1.1)',...
            'pdf: [15,100], t-stat: [1.1,\infty)',...
            'mean \eta: -4.59',...
            'mean m: 1.00',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

end

            
grid on 

plot_eta_results = strcat(path_to_results_folder,'plot_joint_eta_m_FKRB_mixed_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'.png');
saveas(fig1,plot_eta_results) 
hold off

view(3)
plot_eta_results2 = strcat(path_to_results_folder,'plot_joint_eta_m_FKRB_mixed_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_view3.png');
saveas(fig1,plot_eta_results2) 


clear fig1 

%%
% eta
param_results_etafox = strcat('plot_marginaleta_Fox_',Fox_model,'_mixedFEdemoshrt_WTP_v',num2str(subsample),'_inc',num2str(inc),'_sd',num2str(set_seed),'_cty_tau1_grid_v3_all_bs.csv');
param_eta     	  = csvread([ path_to_results_folder param_results_etafox]);
	

tstats2_edge = [0 1.1 max(param_eta(:,6))];
tstats2_color_id = discretize(param_eta(:,6),tstats2_edge);
param_eta_g=[param_eta, tstats2_color_id];

for j=1:2
	group_eta_id(:,j)= tstats2_color_id==j;
end

fig2=figure(2)
hold on
% pdf < 5, tstats < 1.1
h1=scatter(param_eta(group_eta_id(:,1),1),param_eta(group_eta_id(:,1),2),50,'o');
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0.5];

% pdf < 5, tstats >= 1.1
h2=scatter(param_eta(group_eta_id(:,2),1),param_eta(group_eta_id(:,2),2),50,'*');
h2.MarkerFaceColor = [0 0 0];
h2.MarkerEdgeColor = [0 0 0.5];


xlabel('Coefficient on Price: \eta','Fontsize',16)
ylabel('Estimated PDF Weight','Fontsize',16)
ylim( [0,80] );
xlim( [-35,5] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12) 
%a = get(gca,'XTickLabel');
%set(gca,'XTickLabel',a,'FontName','Times','fontsize',18)
 
%line([-5.431913,-5.431913],[0,80],'Color','k' );

if inc==16
 if strcmp(Fox_model,'discrete')   
 line([-5.431913,-5.431913],[0,80],'Color','k' );
 %line([-35,5],[0.7650573,0.7650573],[0,0] ,'Color','k' );
     elseif strcmp(Fox_model,'discretev5')
     line([-5.357534,-5.357534],[0,80],'Color','k' );

     elseif strcmp(Fox_model,'discretev5')
     line([-5.435024,-5.435024],[0,80],'Color','k' );
     
     elseif strcmp(Fox_model,'discretev7')
     line([-5.237717,-5.237717],[0,80],'Color','k' );

     elseif strcmp(Fox_model,'discrete_3tau')    
     line([-5.58872,-5.58872],[0,80],'Color','k' );

     end
 
elseif inc==1
 line([-6.421407,-6.421407],[0,80],'Color','k' );
 %line([-35,5],[0.3539926,0.3539926],[0,0] ,'Color','k' );
elseif inc==2
 line([-6.261386,-6.261386],[0,80],'Color','k' );
 %line([-35,5],[0.4845289,0.4845289],[0,0] ,'Color','k' );
elseif inc==3
 line([-5.604053,-5.604053],[0,80],'Color','k' );
 %line([-35,5],[0.7280716,0.7280716],[0,0] ,'Color','k' );
elseif inc==4
 line([-5.082855,-5.082855],[0,80],'Color','k' );
 %line([-35,5],[0.8185474,0.8185474],[0,0] ,'Color','k' );
elseif inc==5
 line([-5.111113,-5.111113],[0,80],'Color','k' );
 %line([-35,5],[0.9405871,0.9405871],[0,0] ,'Color','k' );
elseif inc==6
 line([-4.592757,-4.592757],[0,80],'Color','k' );
 %line([-35,5],[1.000794,1.000794],[0,0] ,'Color','k' );
end



if inc==16
 if strcmp(Fox_model,'discrete')
 Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.43',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

 elseif strcmp(Fox_model,'discretev5') 
  Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.36',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
 
  elseif strcmp(Fox_model,'discretev6') 
  Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.44',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

 elseif strcmp(Fox_model,'discretev7')
  Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.24',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
 
 elseif strcmp(Fox_model,'discrete_3tau')     
   Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.59',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
 
 end

elseif inc==1
 Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -6.42',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==2
 Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -6.55',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==3
 Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.60',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==4
 Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.08',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==5
 Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -5.11',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==6
Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean \eta: -4.59',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
end
            
grid on 

plot_eta_results = strcat(path_to_results_folder,'plot_marginal_eta_FKRB_mixed_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_v1.png');
saveas(fig2,plot_eta_results) 
hold off

clear fig2

fig3=figure(3)
hold on
% pdf < 5, tstats < 1.1
h1=scatter(param_eta(:,1),param_eta(:,2),50,'o');
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0.5];

for j=1:length(param_eta(:,1))
	line([param_eta(j,1),param_eta(j,1)],[param_eta(j,4),param_eta(j,5)],'Color',[0 0 0.5]);
end


xlabel('Coefficient on Price: \eta','Fontsize',16)
ylabel('Estimated PDF Weight','Fontsize',16)
ylim( [0,80] );
xlim( [-35,5] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12) 
%a = get(gca,'XTickLabel');
%set(gca,'XTickLabel',a,'FontName','Times','fontsize',18)
 

legend('off')
            
grid on 

plot_eta_results = strcat(path_to_results_folder,'plot_marginal_eta_FKRB_mixed_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_v2.png');
saveas(fig3,plot_eta_results) 
hold off

clear fig3


%%
% theta
param_results_thetafox = strcat('plot_marginaltheta_Fox_',Fox_model,'_mixedFEdemoshrt_WTP_v',num2str(subsample),'_inc',num2str(inc),'_sd',num2str(set_seed),'_cty_tau1_grid_v3_all_bs.csv');
param_theta     	  = csvread([ path_to_results_folder param_results_thetafox]);
	

tstats2_edge = [0 1.1 max(param_theta(:,6))];
tstats2_color_id = discretize(param_theta(:,6),tstats2_edge);
param_theta_g=[param_theta, tstats2_color_id];

for j=1:2
	group_theta_id(:,j)= tstats2_color_id==j;
end

fig4=figure(4)
hold on
% pdf < 5, tstats < 1.1
h1=scatter(param_theta(group_theta_id(:,1),1),param_theta(group_theta_id(:,1),2),50,'o');
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0.5];

% pdf < 5, tstats >= 1.1
h2=scatter(param_theta(group_theta_id(:,2),1),param_theta(group_theta_id(:,2),2),50,'*');
h2.MarkerFaceColor = [0 0 0];
h2.MarkerEdgeColor = [0 0 0.5];


xlabel('Misperception: m','Fontsize',16)
ylabel('Estimated PDF Weight','Fontsize',16)
ylim( [0,80] );
xlim( [-3,6] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12) 
%a = get(gca,'XTickLabel');
%set(gca,'XTickLabel',a,'FontName','Times','fontsize',18)
 
%line([0.89,0.89],[0,80],'Color','k' );

if inc==16
  if strcmp(Fox_model,'discrete')  
    %line([-5.431913,-5.431913],[-3,6],[0,0],'Color','k' );
    line([0.7650573,0.7650573],[0,80] ,'Color','k' );
  elseif strcmp(Fox_model,'discretev5') 
     line([0.7911761,0.7911761],[0,80] ,'Color','k' );
  elseif strcmp(Fox_model,'discretev6') 
     line([0.7624894,0.7624894],[0,80] ,'Color','k' );   
  elseif strcmp(Fox_model,'discretev7')
      line([0.7336762,0.7336762],[0,80] ,'Color','k' );
  elseif strcmp(Fox_model,'discrete_3tau')    
      line([0.7778516,0.7778516],[0,80] ,'Color','k' );
  end
 
elseif inc==1
 %line([-6.421407,-6.421407],[-3,6],[0,0],'Color','k' );
 line([0.3539926,0.3539926],[0,80] ,'Color','k' );
elseif inc==2
 %line([-6.261386,-6.261386],[-3,6],[0,0],'Color','k' );
 line([0.4845289,0.4845289],[0,80] ,'Color','k' );
elseif inc==3
 %line([-5.604053,-5.604053],[-3,6],[0,0],'Color','k' );
 line([0.7280716,0.7280716],[0,80] ,'Color','k' );
elseif inc==4
 %line([-5.082855,-5.082855],[-3,6],[0,0],'Color','k' );
 line([0.8185474,0.8185474],[0,80] ,'Color','k' );
elseif inc==5
 %line([-5.111113,-5.111113],[-3,6],[0,0],'Color','k' );
 line([0.9405871,0.9405871],[0,80] ,'Color','k' );
elseif inc==6
 %line([-4.592757,-4.592757],[-3,6],[0,0],'Color','k' );
 line([1.000794,1.000794],[0,80] ,'Color','k' );
end


if inc==16
  if strcmp(Fox_model,'discrete')  
    Lobj=legend('t-stat: [0,1.1)',...
                't-stat: [1.1,\infty)',...
                'mean m: 0.77',...
                'Location','NorthWest');
    EG = findobj(Lobj,'type','text');
    LEG.FontSize = 10;
  elseif strcmp(Fox_model,'discretev5')  
     Lobj=legend('t-stat: [0,1.1)',...
                't-stat: [1.1,\infty)',...
                'mean m: 0.79',...
                'Location','NorthWest');
    EG = findobj(Lobj,'type','text');
    LEG.FontSize = 10;  
  elseif strcmp(Fox_model,'discretev6')  
     Lobj=legend('t-stat: [0,1.1)',...
                't-stat: [1.1,\infty)',...
                'mean m: 0.76',...
                'Location','NorthWest');
    EG = findobj(Lobj,'type','text');
    LEG.FontSize = 10;    
  elseif strcmp(Fox_model,'discretev7')  
    Lobj=legend('t-stat: [0,1.1)',...
                't-stat: [1.1,\infty)',...
                'mean m: 0.73',...
                'Location','NorthWest');
    EG = findobj(Lobj,'type','text');
    LEG.FontSize = 10;
  elseif strcmp(Fox_model,'discretev7')    
    Lobj=legend('t-stat: [0,1.1)',...
                't-stat: [1.1,\infty)',...
                'mean m: 0.78',...
                'Location','NorthWest');
    EG = findobj(Lobj,'type','text');
    LEG.FontSize = 10;
  end    
elseif inc==1
 
Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean m: 0.35',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==2
 
Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean m: 0.48',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==3
 
Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean m: 0.73',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==4
  
Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean m: 0.82',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==5
 
Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean m: 0.94',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
elseif inc==6
 
Lobj=legend('t-stat: [0,1.1)',...
            't-stat: [1.1,\infty)',...
            'mean m: 1.00',...
            'Location','NorthWest');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;
end


            
grid on 

plot_eta_results = strcat(path_to_results_folder,'plot_marginal_theta_FKRB_mixed_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_v1.png');
saveas(fig4,plot_eta_results) 
hold off

clear fig2

fig5=figure(5)
hold on
% pdf < 5, tstats < 1.1
h1=scatter(param_theta(:,1),param_theta(:,2),50,'o');
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0.5];

for j=1:length(param_theta(:,1))
	line([param_theta(j,1),param_theta(j,1)],[param_theta(j,4),param_theta(j,5)],'Color',[0 0 0.5]);
end


xlabel('Misperception: m','Fontsize',16)
ylabel('Estimated PDF Weight','Fontsize',16)
ylim( [0,80] );
xlim( [-3,6] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12) 
%a = get(gca,'XTickLabel');
%set(gca,'XTickLabel',a,'FontName','Times','fontsize',18)
 

legend('off')
            
grid on 

plot_eta_results = strcat(path_to_results_folder,'plot_marginal_theta_FKRB_mixed_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_v2.png');
saveas(fig5,plot_eta_results) 
hold off

clear fig5


%%
%Figure for m comparing FKRB and parametric RCL


rng(matlab_seed)
draws_mixed=500;

	tmp_file1=strcat('/Users/shoude/Dropbox/eegap/EEgap_data_code_heter_SM/Data/Matlab_estimation/choiceset_identifier_trimester_week_zipcode_2008_2012_v11022017_struct_',num2str(subsample),'_',num2str(inc),'_seed_',num2str(set_seed),'.csv');
    id_full=load(tmp_file1); 
	Nmax=size(id_full,1);
    
	if subsample == 11000
		bs_max=10*ceil(Nmax/10000);
    elseif subsample == 44000
		bs_max=10*ceil(Nmax/20000); 
    end

c=1
for j=1:bs_max	
     if subsample == 11000
		index_bs=[(1+(j-1)*(1000)):min(j*(1000),Nmax)];
    elseif subsample == 44000
		index_bs=[(1+(j-1)*(2000)):min(j*(2000),Nmax)];
    end
    
    id=id_full(index_bs,:);
    id(:,1:7)=[];
    choicest_size=sum(id,2);
    %id_discard=find(choicest_size>200 | choicest_size<50);
    id_discard=find(choicest_size>400);
    id(id_discard,:)=[];  
    choicest_size(id_discard)=[]; 
    id=id';
    J=size(id,1);
    
	try	
		result_mixed=strcat('/Users/shoude/Dropbox/eegap/EEgap_data_code_heter_SM/Results_Estimation/demand_est_grID/','mixed_logit_FEdemoshrt_grID_wtp_100iters_',num2str(subsample),'_',num2str(inc),'_',num2str(set_seed),'_matlab_seed_',num2str(matlab_seed),'_nb_draws_',num2str(get_nb_draws),'bs_',num2str(j),'.mat');
        eval(['load ' result_mixed]);
        result_mixed0=theta_est1;
        param_nopt=result_mixed0(1:J-1,1);

        eta(c)=result_mixed0(J)  
        tau(c)=result_mixed0(J+1)
        pi(c)=result_mixed0(J+2)
        theta(c)=result_mixed0(J+3)
        beta(c,:)=result_mixed0(J+4:J+15)
        L_11(c)=result_mixed0(J+16)
        L_21(c)=result_mixed0(J+17)
        L_22(c)=result_mixed0(J+18)
        
        L=[L_11(c),0; L_21(c),L_22(c)];
        mu=[eta(c),theta(c)];
        sigma=L*L';
        R{c} = mvnrnd(mu,sigma,draws_mixed);
	    c=c+1;
	catch
		c=c;	
	end
    
end



for i=1:c-1   
[freq_theta(:,i),edge_theta(:,i)]=histcounts(R{i}(:,2),[min(param_theta(:,1))-0.5;param_theta(:,1)]);
end
id_nz=find(freq_theta>0);
id_pos=freq_theta*0;    
id_pos(id_nz)=1;
weight_theta_mix=100*freq_theta/draws_mixed;
avg_weight_theta_mixed=sum(weight_theta_mix,2)./sum(id_pos,2);

weight_theta_mix(weight_theta_mix==0) = NaN;
weight_theta_std = nanstd(weight_theta_mix');
weight_theta_se = weight_theta_std'./(sum(id_pos,2)).^0.5;

avg_weight_theta_L=avg_weight_theta_mixed-1.96*weight_theta_se;
avg_weight_theta_U=avg_weight_theta_mixed+1.96*weight_theta_se;

param_theta_mixed=[avg_weight_theta_mixed, weight_theta_std', weight_theta_se, avg_weight_theta_L, avg_weight_theta_U,param_theta(:,1)];

[r_z,c_z]=find(param_theta_mixed(:,3)==0);
param_theta_mixed(r_z,:)=[];

fig7=figure(7)
hold on

[r_z,c_z]=find(param_theta(:,3)==0);
param_theta(r_z,:)=[];

h1=scatter(param_theta(:,1),param_theta(:,2),50,'o');
%h1.MarkerFaceColor = [0 0 0.5];
%h1.MarkerEdgeColor = [0 0 0.5];
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0];

h2=scatter(param_theta_mixed(:,6),param_theta_mixed(:,1),50,'*');
%h2.MarkerFaceColor = [1 0 0];
%h2.MarkerEdgeColor = [1 0 0];
h2.MarkerFaceColor = [0 0 0];
h2.MarkerEdgeColor = [0 0 0];
%,'LineWidth',1

for j=1:length(param_theta(:,1))
	%line([param_theta(j,1),param_theta(j,1)],[param_theta(j,4),param_theta(j,5)],'Color',[0 0 0.5],'LineStyle','--');
    %line([param_theta(j,1),param_theta(j,1)],[param_theta(j,4),param_theta(j,5)],'Color',[0 0 0.5]);
 line([param_theta(j,1),param_theta(j,1)],[param_theta(j,4),param_theta(j,5)],'Color',[0 0 0]);

end


for j=1:length(param_theta_mixed(:,6))
    %if avg_weight_L(j)~=NaN;
        %line([param_theta_mixed(j,6),param_theta_mixed(j,6)],[param_theta_mixed(j,4),param_theta_mixed(j,5)],'Color',[1 0 0]);
    line([param_theta_mixed(j,6),param_theta_mixed(j,6)],[param_theta_mixed(j,4),param_theta_mixed(j,5)],'Color',[0 0 0]);        %end
end    


xlabel('Misperception: m','Fontsize',16)
ylabel('Estimated PDF Weight','Fontsize',16)
ylim( [0,80] );
xlim( [-3,6] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12) 
%a = get(gca,'XTickLabel');
%set(gca,'XTickLabel',a,'FontName','Times','fontsize',18)
 
if inc==16

    Lobj=legend('FKRB, E[m]=0.77',...
            'RCL, E[m]=0.64');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

end            
grid on 

plot_theta_w_mixed_results = strcat(path_to_results_folder,'plot_marginal_theta_FKRB_mixed_vs_mixedonly_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_v2.eps');
saveas(fig7,plot_theta_w_mixed_results) 

clear fig7


%Marginal for eta

for i=1:c-1   
[freq_eta(:,i),edge_eta(:,i)]=histcounts(R{i}(:,1),[min(param_eta(:,1))-0.5;param_eta(:,1)]);
end
id_nz=find(freq_eta>0);
id_pos=freq_eta*0;    
id_pos(id_nz)=1;
weight_eta_mix=100*freq_eta/draws_mixed;
avg_weight_eta_mixed=sum(weight_eta_mix,2)./sum(id_pos,2);

weight_eta_mix(weight_eta_mix==0) = NaN;
weight_eta_std = nanstd(weight_eta_mix');
weight_eta_se = weight_eta_std'./(sum(id_pos,2)).^0.5;

avg_weight_eta_L=avg_weight_eta_mixed-1.96*weight_eta_se;
avg_weight_eta_U=avg_weight_eta_mixed+1.96*weight_eta_se;

param_eta_mixed=[avg_weight_eta_mixed, weight_eta_std', weight_eta_se, avg_weight_eta_L, avg_weight_eta_U,param_eta(:,1)];

[r_z,c_z]=find(param_eta_mixed(:,3)==0);
param_eta_mixed(r_z,:)=[];

fig8=figure(8)
hold on

[r_z,c_z]=find(param_eta(:,3)==0);
param_eta(r_z,:)=[];

h1=scatter(param_eta(:,1),param_eta(:,2),50,'o');
%h1.MarkerFaceColor = [0 0 0.5];
%h1.MarkerEdgeColor = [0 0 0.5];
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0];

h2=scatter(param_eta_mixed(:,6),param_eta_mixed(:,1),50,'*');
%h2.MarkerFaceColor = [1 0 0];
%h2.MarkerEdgeColor = [1 0 0];
h2.MarkerFaceColor = [0 0 0];
h2.MarkerEdgeColor = [0 0 0];
%,'LineWidth',1

for j=1:length(param_eta(:,1))
	%line([param_theta(j,1),param_theta(j,1)],[param_theta(j,4),param_theta(j,5)],'Color',[0 0 0.5],'LineStyle','--');
    line([param_eta(j,1),param_eta(j,1)],[param_eta(j,4),param_eta(j,5)],'Color',[0 0 0]);
end


for j=1:length(param_eta_mixed(:,6))
    %if avg_weight_L(j)~=NaN;
        line([param_eta_mixed(j,6),param_eta_mixed(j,6)],[param_eta_mixed(j,4),param_eta_mixed(j,5)],'Color',[0 0 0]);
    %end
end    


xlabel('Coefficient on Price: \eta','Fontsize',16)
ylabel('Estimated PDF Weight','Fontsize',16)
ylim( [0,60] );
xlim( [-35,5] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12) 
%a = get(gca,'XTickLabel');
%set(gca,'XTickLabel',a,'FontName','Times','fontsize',18)
 
if inc==16

    Lobj=legend('FKRB, E[\eta]=-5.43',...
            'RCL, E[\eta]=-4.98');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

end            
grid on 

plot_theta_w_mixed_results = strcat(path_to_results_folder,'plot_marginal_eta_FKRB_mixed_vs_mixedonly_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_v2.eps');
saveas(fig8,plot_theta_w_mixed_results) 

clear fig8



%Joint for eta and theta
edges_fox{1}=[-35,-30,-25,-18,-16,-14,-10,-7,-5,-4.5,-4,-3.5,-3.25,-3,-2.5,-2,-1,-0.1,1,2,3];
edges_fox{2}=[-3.5,-2.5,-2,-1.5,-1.25,-1,-0.75,-0.5,-0.25,-0.1,0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1,1.1,1.25,1.5,2,2.5,3,3.5,4,5];

for i=1:c-1           
        freq_mixed=hist3(R{i},'Edges',edges_fox);
        weight_mixed(:,i)=100*freq_mixed(:)/draws_mixed;
end
id_nz = find(weight_mixed>0);
id_pos = weight_mixed*0;    
id_pos(id_nz) = 1;
avg_weight_mixed = sum(weight_mixed,2)./sum(id_pos,2);

weight_mix(weight_mixed==0) = NaN;
weight_std = nanstd(weight_mixed');
weight_se = weight_std'./(sum(id_pos,2)).^0.5;

avg_weight_L=avg_weight_mixed-1.96*weight_se;
avg_weight_U=avg_weight_mixed+1.96*weight_se;

% meshgrid...
[x_eta,y_theta]=meshgrid(edges_fox{2},edges_fox{1});
%cc=[aa(:),bb(:)];
param_mixed=[avg_weight_mixed, weight_std', weight_se, avg_weight_L, avg_weight_U,x_eta(:),y_theta(:)];

[r_z,c_z]=find(param_mixed(:,3)==0);
param_mixed(r_z,:)=[];

fig9=figure(9)
hold on

[r_z,c_z]=find(param_eta(:,3)==0);
param_eta(r_z,:)=[];

h1=scatter3(param_mixed(:,6),param_mixed(:,7),param_mixed(:,1),50,'o');
h1.MarkerFaceColor = [0 0 0];
h1.MarkerEdgeColor = [0 0 0.5];

xlabel('Coefficient on Price: \eta','Fontsize',16)
ylabel('Misperception: m','Fontsize',16)
zlabel('Estimated PDF Weight','Fontsize',16)
ylim( [-3,6] );
xlim( [-35,5] );
%xt = get(gca, 'XTick');
%set(xt, 'fontsize', 16)
set(gca,'fontsize',12)   
 
if inc==16

    Lobj=legend('FKRB, E[\eta]=-5.43',...
            'RCL, E[\eta]=-4.98');
EG = findobj(Lobj,'type','text');
LEG.FontSize = 10;

end            
grid on 

plot_eta_theta_w_mixed_results = strcat(path_to_results_folder,'plot_joint_eta_theta_FKRB_mixed_vs_mixedonly_',Fox_model,num2str(subsample),'_inc_',num2str(inc),'_seed_',num2str(set_seed),'_v2.png');
saveas(fig9,plot_eta_theta_w_mixed_results) 

clear fig9


end